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Abstract 

We present an efficient Monte Carlo method to study the properties of the 
infinitely repulsive Hubbard model with few holes at finite nonzero temper- 
atures. Using this method, magnetic and specific heat results are presented 
for a 10 x 10 lattice with two holes. Our results show that the magnetization 
increases with decreasing temperature for such a finite system. 

I. INTRODUCTION 

The magnetic behaviour of the Hubbard model in the U = 00 limit is an old problem. 
For this case a well-known exact theorem is due to Nagaoka (1966) and Thouless (1965). 
It provides explicit, but highly idealized conditions under which ferromagnetism is stable in 
the Hubbard model. It proves that for infinitely large repulsive Hubbard U the macroscopic 
degeneracy of the ground state at half filling (when the number of electrons N is exactly 
equal with the number of lattice sites, hereafter denoted by N\) is lifted by a single hole, 
i.e., when N = N\ — 1. In this case the saturated ferromagnetic ground state is stable 
for any value of the hopping amplitude t on square, simple cubic and bcc lattices, and for 
t < on triangle, fee and hep lattices. A very elegant and fully general proof was given 
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for this result by Tasaki (1989) based on the Perron-Frobenius theorem. For the Nagaoka 
mechanism to work, the lattice needs to contain loops along which the hole can move. Once 
the hole moves, the maximal overlap between the initial and the final state clearly occurs in 
a ferromagnetic configuration. The problems are that Nagaoka's proof does not even extend 
to two holes, that a single hole is thermodynamically irrelevant, and that the limit of U — oo 
is highly unrealistic. Moreover, several theoretical results state the instability of Nagaoka's 
ferromagnetism for two holes, see e.g. Doucot and Wen (1989), Siito (1991 a), Toth (1991). 

In a recent published paper (Gurin and Gulacsi, 2000) a formalism was presented which 
helps to deduce the partition function and some physical quantities of the system based on 
a loop summation (i.e. path integral on lattice) technique. The presented formulas have 
the advantage that the trace over the spin degrees of freedom has been taken exactly into 
account. Resembling technique was used by Siito (1991 b) to analyze the magnetic behaviour 
of the system in thermodynamic limit. However, his result is valid in thermodynamic limit, 
only in the case N h / N A — > where N h = N A — N and not for finite hole concentration 
and T — > 0. In this paper we investigate the behaviour of a small system with finite hole 
concentration, using Monte Carlo method at T ^ 0. Starting from the aim to enhance the 
study of magnetic properties in the presence of a small hole concentration even at finite 
nonzero temperatures, the procedure is described in detail. The studied system is a 10 x 10 
lattice with periodic boundary conditions containing two holes. The results show an increase 
of the magnetisation with the decrease of the temperature. 

The remaining part of the paper is organized as follows: in Sec. II. we describe the 
representation of the Hilbert space in which the calculation is done, Sec. III. presents the 
partition function and the analyzed physical quantities, and Sec. IV. contains the Monte 
Carlo results obtained, together with a discussion of the deduced results. 

II. REPRESENTATION OF THE HILBERT SPACE 

Our Hamiltonian describes the one-band Hubbard model in U — oo limit: 
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A.. = y: Uj H% = £ty p ( E ^ )p, (i) 

where Uj = — t are hopping matrix elements for nearest neighbour sites (otherwise Uj = 0). 
The operators c| jCr ( Cj jCr ) creates (annihilates) an electron with spin a at lattice site i, and 
hi = So- ^i,o- represents the particle number operator. The double occupancy is projected 
out by P = J2i ( 1 - ). 

In our study we need operations with permutations therefore we fix below some notations 
concerning them. We treat a permutation V of degree n as a bijective function of the set of 
the first n natural numbers into themselves. The product of the permutations is the standard 
function composition: ( V Q )( i ) = V{ Q( i ) ). The symmetric group of degree n formed 
by all the permutations of degree n is denoted by S n . The complex Hilbert space generated 
by the elements of S n as an abstract orthonormalized basis endowed by the product defined 
by the convolution is called the group algebra of the symmetric group. The parity of a 
permutation is | V |. 

We study the one band U — oo Hubbard model. In such a tight-bindind description the 
electrons can live only on lattice sites in two orthogonal states with spin up and down. Due 
to the kinetic energy term, an electron hops from one site to an other (while its spin state 
remains unchanged), changing however the locations of the empty and occupied sites. Here- 
after, fixed positions of the occupied sites we call charge configuration which does not give 
account of the electronic spin states, only their places. Iterating the effect of the hopping 
term on a state several times, hopping around the lattice, the electrons produce different 
charge configurations, particulary they can produce the original configuration again as well. 
However, the spin up and spin down electrons are commingled, therefore the spin configu- 
ration can be different than in the original state. In our descriptiom a spin configuration 
is given by the position of sites occupied by spin up and separately spin down electrons, 
respectively. In this way, fixing the charge configuration and the third component of the 
total spin S z , the spin configuration can be described by an ordered set of the up and down 
spin states: (cri)f =1 in which <7; G { |, | } indicates the spin state of the electron situated 
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in the i-th occupied lattice site. Certainly, at first, the occupied sites have to be labelled 
from 1 to N. In this way, every spin configuration can be obtained by a permutation from 
(<7°(0)£i = (T, T,...,T, I, I, •••,!), where a°(i) = T if i = 1, . . . , iV T = N/2 + S* and 
a°(i) = jif i = Ni + 1,...,N. 

Based on the above mentioned observations, the Hilbert space of the system and a basis 
of it can be expressed in terms of permutations. The charge configuration can be described 
by a function K : { 1 , ... , N A } -> { 1 , ... , N\ } for which 11(1), ... , K(N) are the 
occupied lattice sites (taken into account an arbitrary but fixed numbering of sites), and 
1Z(N + 1), ... , TZ(N\) are the positions of holes. Actually 1Z G Sn a , but hereafter we 
denote the charge configuration by Tl for which we require that 72.(1), ... , TZ(N) have to 
be increasingly ordered. Moreover, we do not distingwish permutations which differ only in 
the order of the positions of the holes 1Z(N + 1), ... , 1Z(N\), everyone is denoted by 7Z, and 
the further definitions and expressions do not depend on this undetermination. The Hilbert 
space generated by all the different charge configurations as an orthonormalized basis is 
denoted by 7i c . 

To label the occupied sites from 1 to N we use the function V\r l . After it, a spin 
configuration ( o~i )f =1 can be described by a permutation V G Sjy which rearrange the basic 
order in (cr°(i) )fL 1 in such a way that a°(V(i) ) = Oi. This means that in a state described 
by a charge configuration TZ and a spin configuration V, the electron situated in the lattice 
site i is in the spin state <x°( V'R,^ 1 ( i)). More precisely, there is no unique connection 
between the spin configurations and permutations. A permutation uniquely determines a 
spin configuration, but the reverse is not true. Let Q G Sjy be a permutation which 
rearranges the elements of the subset (a°)~ 1 ( j ) = { 1, . . . , }, and rearranges the 
elements of (o-°) _1 ( | ) = { + 1, . . . , N }, but does not interchange the elements of these 
two subsets. It is clear for arbitrary permutation V, that V and QV describe the same spin 
configuration, because <t°( V(i)) = cr°( QV(i) ) for all i. Mathematically this means 
the following. The above defined permutations Q form a subgroup Kg* C Sn, and a right 
coset V := { QV : Q G Ks* } is in one-to-one correspondence with a spin configuration. 
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Moreover, the Hilbert space generated by these right cosets as an orthonormalized basis, is 
isomorphic with a right ideal of the group algebra, which will be denoted by 7Yf z . On this 
line, an arbitrary state with given S z , can be written as ip c <g> <£>f z where (p c G Ti. c and 
y?f z G 7if z . Furthermore, {|72) } is a basis of H c , and {\P) } is a basis of 7if z , respectively 
The explicit expression of the above defined basis with the usual creation operators is 

N 

I n > ® | v ) = n 4(0^(0) 1 > > ( 2 ) 
i=i 

where P is an arbitrary permutation of the right coset V. Now it can be seen that the 
requirement that 72.(1), 7Z(N) have to be increasingly ordered is needed to fix the order of 
the creation operators in Eq.(Q). Let C a permutation which rearranges the first N numbers 
and leaves unchanged the last Nh ones. Using 1Z' = 1ZC instead of 72 in Eq.(|]), we get the 
same operators but in different order, supposing that it is also 72' which is used to label the 
occupied sites from 1 to N. This rearrange of the creation operators involves a sign (— 1)' C L 
To avoid this sign ambiguity, we use 72 in Eq.(fJ), for which the order of the operators is well 
determined due to the requirement for the order of 72.(1), ...,TZ(N). 

III. PARTITION FUNCTION AND PHYSICAL QUANTITIES 

Now we study the changing of a state vector from Eq.(U) during the motion of the 
particles. Let's study the hopping of particles step-by-step. Suppose, that the starting 
point is 72 (g> V and in the first step an electron hops from site i\ to site ji, in the next step 
an other electron (or the same) from 12 to j'2 and so on. The path in the charge configuration 
space 7 = ( 72 , 72i , . . . , 72;„i , 72/) is uniquely determined by the sequence of nearest 
neighbour pairs ( %\ , ji ) , ... , (ii , ji)- In the nth step a hole and an electron interchange 
their positions. This can be described by the transposition 7 >in ' Jn . Now we suppose, that 
after I steps the system gets back to the original charge configuration. However the positions 
of holes are the same as in the original situation, the moving rearranged the electrons. The 
electron situated originally at site i moved to site j. The new situation can be described by a 
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permutation TZ^TZ where TZ 1 = V H,jl . . . V' ll,jl G Sn a which takes % to j for all % which label 
occupied sites and permutes somehow the labels of empty sites. Because the spin state of 
the electron remains unchanged during the hopping, the electron situated after the motion 
in site j = TZ^i) is in the spin state a] = a t = <7 (Pft -1 (i)) = ^(VTZ^TZ^j)). The 
charge configuration is the same as the starting TZ therefore to recognize the final spin state 
V' we have to write the final spin state of the electon situated in site j as cr°( V' 1Z~ l (j ) ). 
As a consequence, we have V' — V TZ^ 1 1Z 1 TZ = V V~ l where P 7 = TZ^ 1 1Z~ l TZ is an 
element of Sn a but it rearranges the first N numbers among themself and leaves unchanged 
the last N h ones, therefore it can be considered as an element of SV Furthermore the final 
charge configuration is again TZ = ( TZ 1 TZ ) P 7 instead of TZ y TZ therefore we have to take 
into account a sign (— l)^ — as we mentioned at the end of the previous Section — to 
write the final state as TZ <g> (-l)l'M V V' 1 

The mapping V — > (— l)'^ 7 ' "P V~ 1 is defined on the basis of Tl s * ' . The influence caused 
by a hopping of holes along a loop 7 is described by a linear operator which is the linear 
extension of the mapping V -»• (—1)1^1 V V- 1 to the whole Tif . It is actually a linear 
representation of the symmetric group Sn on 7if z (Gurin and Gulacsi, 2000) hereafter 
denoted by T SZ [V^]. 

Based on the above considerations, we can give expressions for the partition function as 
follows 

Z = Tr e~— = Tr « sz e ~~ , (3) 

s z 

and the expectational value of the square of the total spin 

(S 2 ) = ^Tr(5 2 e-#) = 1^3(^) 2 Tr^e-# (4) 

in which the sum over the spin degrees of freedom is expressed analitically in terms of the 
characters of the symmetric group. From the above expressions it can be seen that what we 
need for calculations is in fact 

z=o 1 ■ KK 1 > a 7en^(i) v 
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where f2^( I ) denotes the set of the Z-length loops (closed pathes) with starting point TZ 
in the space of charge configurations, and x is the character of the representation T sz . 
Every character is constant on an arbitrary conjugate class C . Therefore, the above sum 
over fl^( I ) has (I) identical members, where (I) represents the number of the 
loops with length I and starting point TZ for which V 1 G C. Now we use (see Gurin and 
Gulacsi, 2000) 

N/2 

E x s * (c) = (-i)ioi 2 ££x<\ (6) 

S*=-N/2 

furthermore 

E (^)V^) = ("1) |C| 2^ ( jE^Cr • (7) 

S z =-N/2 V «=1 / 



Here | C | is the parity of permutations from the conjugate class C and ( Ci , Ci , . . . , Cjv ) 
describes their cycle structure, i.e. these permutations contain Ci cycles with length i. 
Finally, inserting Eqs.(|^J^,0) into Eqs.(|3],f|) we obtain 




To determine the coefficients (I) we should follow the following procedure. Starting 
from a charge configuration TZ, holes moves / steps in such a way that finally the sysyem 
goes back to the original hole configuration (this is a loop with length /). In this process, 
an arbitrary hole can be moved in every step, the order of the steps being relevant. A step 
means that we interchange a hole with an electron from its neighbourhood. Counting the 
number of Z-loops for which the cycle structure of the permutation of electrons after the Zth 
step is C, we get the number (I). 

7 



IV. MONTE CARLO RESULTS AND DISCUSSIONS 



The coefficients N£>(1) can be determined by Monte Carlo method. By random sampling 
we can obtain the percentage of loops with a given length for which V 7 G C . The Eqs.(||) 
and (||) have the advantage that the trace over the spin degrees of freedom has been taken 
already. Hence we need only to sample the loops (i.e. "world lines") of holes instead of the 
loops of the electrons. This is advantageous since due to the presence of different possible 
electron spin configurations, the electron loops have 2 N times more starting points than the 
hole loops. 

We studied by this method a 2D 10 x 10 lattice with periodic boundary conditions in the 
presence of two holes. We took into account every possible starting hole (charge) configu- 
ration. Because of the periodic boundary conditions, translational, rotational and reflection 
symmetries of the lattice, there are 20 essentially different charge configurations instead of 
( 1 2°) • Because of the fact that the contributions of the 2 98 different spin configurations cor- 
responding to a fixed charge configuration are taken into account analitycally by Eqs.(||||), 
we have taken into account every possible strating point of loops. Furthermore, starting 
from a charge configuration, we took into consideration every path up to the length I = 14. 
Thus, we counted exactly every loop up to / < 14, therefore the coefficients of our high 
temperature expansion are exact up to 14th order in (3 / \ t \ . Then, we continued every path 
with random directions up to / = 100. 

The results can be seen in the Fig. 1. The continuous lines show the result of the 
high temperature expansion with exact coefficients up to 14th order. Taking into account 
further coefficients determined by Monte Carlo method up to 50th order we get the short- 
dash lines, and up to 100th order we get the long-dash lines. The expectation value of the 
square of the total spin goes to 3 iV / 4 = 73.5 (and not to zero) as T goes to infinity. 
It is understandable, because at very high temperature the likelyhood of every possible 
states become equally 1 / Z. (The spectrum of the Hamiltonian is bounded.) Therefore, the 
expectation value of an arbitrary quantity goes to its simple algebraic average over the all 
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states. This gives the result mentioned above for ( S 2 )). This means in fact paramagnetic 
behaviour, because a smaller total spin value has greater thermodynamic weigh, i.e. the 
dimension of the eigensubspace of S 2 is greater, therefore the S = subspace is the greatest 
one, so this is the most probable value of the spin. Because of this reason, the total spin 
(but not ( S 2 )) per particle is proportional to 1 / which goes to zero in thermodynamical 
limit. 

Because of the inaccuracy and the complete neglection of the higher order coefficients in 
our expressions, the numerical error increases as the temperature decreases. Therefore we do 
not plot the T — > temperature region in Fig. 1. Nevertheless, in the presented region the 
results are correct since the higher order terms do not improve the lower order contributions 
in this domain. At these temperatures however, due to the Monte Carlo method itself, small 
numerical errors are present of about 1%. 

The obtained results show that the spontaneous magnetization increases as the tempera- 
ture decreases. This tendency becomes stronger when we use more accurate approximation. 
At the same temperature range when ( S 2 ) starts to increase, the specific heat has a max- 
imum. Hovewer this peak in the specific heat seems to be significant — in higher order 
approximations the curves have to break down as well since the specific heat has to go to 
zero in T — > limit, — but the increase in ( S 2 ) is small, nevertheless very fast. 

Related to the two and few hole cases in the literature it have been published (Doucot and 
Wen 1989, Suto 1991 a, Toth 1991) that the ground state can not have maximal spin. This 
statement is based on variational calculations. There exists trial state with spin S = S max — 1 
which have lower energy than the fully polarized state with spin S max . However, up to now 
nobody can construct variational state with spin S < S max — 2 which would have lower 
energy, thus these results do not exclude the ferromagnetic behaviour, only the highest spin 
value. Simultaneously our results suggest that the very low spin values can be exclude as 
well, and the ground state is not a singlet. As a conclusion we may state that the system 
behave at low temperature as a not fully saturated ferro magnet. This is consistent with 
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exact diagonalization results for small clusters (Arrachea 2000). 
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FIGURES 

FIG. 1. (a): square of the total spin and (b): specific heat as the function of the temperature 
measured in \ t\ units. The continuous (lower) lines show the result of a 14th order high temperature 
expansion, the short-dash lines (in the middle) show the result of an 50th order high temperature 
expansion, the long-dash (upper) lines of a 100th order expansion based on coefficients determined 
by Monte carlo method. 
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